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Abstract. We propose and test a new method for generating canonical sequences 
for analysis by the Coherent Anomaly Method (CAM) from non-mean-ficld 
approximations. By intentionally underestimating the rate of convergence of exact- 
diagonalization values for the mass or energy gaps of finite systems, we form families 
of sequences of gap estimates. The gap estimates cross zero with generically nonzero 
linear terms in their Taylor expansions, so that v= 1 for each member of these sequences 
of estimates. Thus, the CAM can be used to determine v. Our freedom in deciding 
exactly how to underestimate the convergence allows us to choose the sequence that 
displays the clearest coherent anomaly. We demonstrate this approach on the two- 
dimensional ferromagnetic Ising model, for which v=\. We also use it on the three- 
dimensional ferromagnetic Ising model, finding v w 0.629, in good agreement with 
other estimates. Finally we apply it to an antiferromagnetic spin-1 Heisenberg chain, 
finding vk, 0.987 at the phase transition between the Haldane phase and the dimerized 
phase, in agreement with the field-theoretic prediction v = 1. Although the specific 
systems used to test the extrapolation- CAM procedure involve finite system sizes, the 
method could be applied to other finite approximations, such as systematic variational 
approximations. 



PACS numbers: 64.60.Fr, 05.50.+q, 05.30.-d 
1. Introduction 

On approaching a critical point, some quantity diverges in the thermodynamic limit 
with a characteristic critical exponent. The Coherent Anomaly Method (CAM) has 
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proven quite successful in determining critical exponents from certain sequences of 
approximations [1, 2, 3]. The CAM requires a systematic, or canonical, sequences of 
approximations, all of which yield identical, known critical exponents. The prototypical 
example is a sequence of mean-field approximations in which successively larger clusters 
allow more and more fluctuations to be properly taken into account [4]; all of the 
critical exponents in this example assume their "classical" values. The CAM uses the 
known critical exponents of the approximate systems and the behaviour of the critical 
amplitudes as the quality of approximation is improved to determine the true critical 
exponents of the original system being approximated. The purpose of this paper is 
to show how sequences of approximations in which the critical points are ill defined 
and divergences do not occur can be used to construct canonical sequences through 
extrapolation. 

The basic idea is as follows. Suppose that in the thermodynamic limit a nonnegative 
quantity £((3) diverges as a system parameter f3 approaches a critical value f3* as 

m ~ i& - p)- v , a) 

where v is a critical exponent of unknown value. The reciprocal of £(/3) then converges 
to zero as 

m^wr 1 ~{&-er- ( 2 ) 

Suppose further that we have a sequence of monotonically decreasing approximations 
{Si(/3)} such that 

W) > 5 i+1 (f3) > <U/3) = A(/3) V(i,/3) . (3) 

By taking two or more consecutive values of 5i(ft), we make an extrapolation A{j}(/3) 
which intentionally underestimates the rate of convergence in {5i([3)}. Clearly, these 
extrapolations will be negative at P*, and they generically cross zero as 



A W (/3) 



dA 



{<} 



dp 



P = P{i},c 



(%,c - P) (4) 



for some {/%}, c }- For any reasonable extrapolation procedure (e.g. power-law 
extrapolation or exponential extrapolation), lim{j}_ >0O A^y(P) = A(/3). The CAM 
hypothesis, which can be justified on the basis of an envelope argument [4], is that 



dA 



- (p: - /%,cr 1 , (5) 

P = /%},c 



dp 

so that, in analogy with (2), 

A {i} (p) ~ (p* - p {i}>c y-\p {i}>c - p) . ( 6 ) 
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This provides us with a convenient means of measuring v. According to (5), a plot of 

dA w 

P = /%} 



Y {l} = In -/%, 



d/5 



(7) 



versus 



should (in the limit X^—oo) be a straight line with slope v—\. Furthermore, since all 
that is required of the extrapolation is that it must underestimate the convergence, we 
are free to choose an extrapolation which leads to a particularly clear coherent anomaly. 

The organization of the remainder of this paper is as follows. In section 2 we take as 
A(/3) the mass gap, i.e., the reciprocal of the correlation length, as a function of inverse 
temperature in the square-lattice Ising ferro magnet. This model has the advantage 
that the correlation length is known analytically [5]. In section 3 we analytically study 
the asymptotic behaviour of the estimated critical exponent if Si is given by a finite-size 
scaling function. In section 4 we study the critical behaviour of the mass gap in the cubic- 
lattice Ising ferromagnet, and demonstrate that the method works even when v^fi 1. In 
section 5 we study the critical behaviour of the energy gap in an antiferromagnetic spin- 
1 Heisenberg chain with bilinear and biquadratic interactions at the phase transition 
between the Haldane and dimerized phases. In section 6 we summarize and discuss 
possible extensions of this work. 

Note that although all of the initial approximations 5i used in this study come from 
systems that are finite in at least one dimension, this is only to provide convenient 
examples. The method itself does not restrict us to such systems. 



2. The Square-Lattice Ising Ferromagnet 

As the first example, we use the extrapolation-CAM method to determine v for the 
classical Ising ferromagnet with Hamiltonian 

n = - E > ( 9 ) 

where Sj = ±l. The sum Yuj) runs over a U nearest-neighbour pairs on a periodic square 
lattice which is of length L in the y-direction and of infinite length in the x-direction. 
The unit of length is the lattice constant. 

This model has some very advantageous properties: the mass gaps {5l((3)} can 
be calculated analytically for systems of arbitrary finite width L, and the mass gap 
A(/3) for the thermodynamic limit of the model can also be calculated analytically [5]. 
Consequently, we can make a detailed comparison of the extrapolation-CAM estimates 
for the critical exponent v with its rigourously known value, v — \. 
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For this model the variable (3 is the reciprocal of the dimensionless temperature. 
Since it is known that 5 L ~L~ l for sufficiently large systems at the critical point [6, 7, 8], 
we extrapolate by solving 



for A i L /(/3) at fixed (5 for B<1. The result of such an extrapolation is shown in figure 1. 

Figure 2 shows Y vs. X, defined by (7) and (8), to be a curve with a small slope. 
The slope tends to zero as —X becomes large — that is, when Pl,l', c becomes a very 
good approximation to (3*. We can estimate v from the slope of the line connecting two 
adjacent points (Xi, Yi) and (X 2 , Y 2 ). This estimate will depend not only on the quality 
of the initial approximations 5l(/3), but also on the parameter B. We can constrain B 
by taking a third point (X 3 , Y 3 ) and demanding that the three points be colinear, so 
that the value of vsj^ will be unambiguous. Figure 3 shows this for values of v based 
on systems with L = 4, 9, 16, and 25. The resulting estimate is =0.987406, in good 
agreement with the exact value v — \. As may be expected, increasing the sizes of the 
four systems needed to form the estimate of v increases the accuracy of V{l}- Figure 4 
shows the estimated value of v vs. B for L = j 2 , (j + l) 2 , (j + 2) 2 , and (j + 3) 2 , where 
j £ {2, 3,4,..., 20}. The estimate for j = 20 is v = 1.000 004. Further increase of j 
actually causes the accuracy of the estimate to become worse due to the increasingly 
large sums required to calculate 5l(/3) and d^/d/? and the finite (eight-byte) numerical 
precision of our programs. The convergence of z/r^} depends on how the four system 
sizes are chosen, and can be both complicated and nonmonotonic. 

3. Relation to Finite-Size Scaling 

Some insight into the dependence of z/{l} on the system sizes can be gained by assuming 
that 5 L satisfies the scaling equation [8] 



S L (P) =A L , L ,([3) + A([3)L 
6 l ,(P)=A LiL ,(P) + A(P)L 



-B 



l-B 



(10) 



6 L (I3) = L"Q(x,y) , 



(11) 



where 




(12) 



and 



(13) 



This assumption is well justified for the systems actually analyzed in this paper, since 
they all are finite in at least one dimension. In order to be consistent with (2), we must 
have ujjQ = —u, and = y T [8]. The variable ( is the correction-to-scaling amplitude 



5 



and is the correction-to-scaling exponent, where the correction to scaling is assumed 
to arise from the leading irrelevant field [9]. 
For large L' = L+1, (10) yields 



L~ B - L'- B 



dQ 



(B + u)Q(x,y) + ex^(x,y) 



+ <Pyf-(x,y) + 0(L- 1 ) 
oy 



(14) 



For the remainder of this section, we consider L to be a single continuous variable and 
drop the second index V . The 0(L _1 ) term in (14) comes from the truncation error 
in a difference formula approximation of a derivative. This term cannot be neglected 
if < —1, unless we generalize the extrapolation procedure to allow for a sufficiently 
high-order difference formula [10]; for the remainder of this section we assume that this 
is done if necessary and neglect the truncation error. 

For the moment, we set C = and perform the extrapolation-CAM procedure in the 
absence of corrections to scaling. 

It is convenient at this point to use a slightly different definition for Y, 



y = in -p: 



dA, 



Using (14), we find 

Y = In 



x 1 



dp 



(15) 



dA, 



dx 



+ [lnx-X]-\nB + ln—(x,B). 

We find Pl,c from the condition A L (p LjC ) — 0, which implies 

dO 

F(x, B) = (B + u)Q(x, 0) + 6x-^-(x, 0) = . 

ox 



(16) 



(17) 



The estimate for v is given by calculating dYjdX while holding both F(x,B) and B 
constant, which is accomplished by holding x constant: 

dF 



dX 



Finally, the "local straightness" constraint 
d 2 Y 



dX 2 







18) 



(19) 
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is simply the identity = 0. Thus without corrections to scaling, the v L = v but x and 
B are undetermined. Note, however, that if x = 0, B = —lu> 0. 

In order to study the effects of corrections to scaling, we could expand (14) and 
the left-hand side of (19) (which is proportional to () in x, y, and e = B+ou; then the 
requirements Al(Pl,c) = and d 2 F/dX 2 = yield a pair of equations of the form 

M n x + M 12 e = v lV ^ 

M 2 \x + M 22 e = v 2 y 

where M and v are constants. The solution of (20), 

= M~Vy, (21) 




shows that both x and e are proportional to y. In the same way, we can expand the 
left-hand side of (18) in terms of x, y, and e. The result then is that 

B ~ -co + 0(L*) (22) 

and 

v L ~v + 0(L*). (23) 

Privman and Fisher have shown that the asymptotic convergence for finite-size scaling 
renormalisation techniques is also of the form (22) [11]; their discussion of the difficulties 
in actually observing this asymptotic behaviour should be relevant in the present case as 
well. For the two-dimensional Ising model, uj = —1, so we should expect lim^oo B — l, 
as indeed seems plausible from figure 4. 



4. The Cubic-Lattice Ising Ferromagnet 

In three dimensions, the Ising model defined by (9) has not been solved analytically, 
but it has been the subject of a large amount of numerical study. We use a Monte Carlo 
Renormalisation Group estimate for the critical point of the cubic-lattice Ising model, 
/3* = 0.221 652(4) [12]. Recent estimates of v include v = 0.642(2) [12], v « 0.631 [13], 
and i/w 0.646 [14]. 

As a result of some earlier studies [15, 16], we have transfer-matrix results already 
available for cubic-lattice Ising ferromagnets with periodic boundary conditions and 
square LxL cross-sections, where L G {2,3,4,5}. These system sizes are obviously 
quite small, and are not really competitive with some current methods, such as the 
Transfer-Matrix Monte Carlo method [17]. 

The method here is basically the same as in the previous section, although we are 
much more restricted both in the system sizes available and, for L = 5, the number of 
data available. Since references [15, 16] dealt with phenomena for (3>{3*, we have only 
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a few points S 5 (P) with (3<{3* [see figure 5]. We use cubic splines to evaluate both S L ((3) 
and the first d5 L /d(3 as continuous functions. Except at the low-/? end of the L = 5 data, 
we always use "natural" boundary conditions on our cubic splines, i.e. specifying that 
d 2 <5i/d/3 2 = at the end points of our data. For the low-/? end of the £5 spline, we use 
both natural boundary conditions and "clamped" boundary conditions by extrapolating 
6.65 /dp from the smaller system sizes. There is little difference in 5 5 itself for these two 
splines. 

We extrapolate A LjL i using (10) and use (5), (7), and (8) to estimate v. As figure 6 
shows, the splines yield unreliable estimates near the end of the £5 data. However, as 
long as /?{4,5}, c is reasonably large, the estimates of v do not depend on the boundary 
conditions used for the spline and seem more reliable. Choosing five such points and 
using Aitken's A 2 method [18] t to accelerate the convergence, we extrapolate to 
[z/(2,3;3,4) — z/(3, 4; 4, 5)] _1 — > 00, yielding v = 0.629. Given the very small systems 
used in this estimate, this is in good agreement with other recent estimates of v. 

5. The Spin-1 Antiferromagnetic Heisenberg Chain 

The Heisenberg chain we study is defined by the quantum Hamiltonian 



where Si is the quantum spin-1 operator for the spin at site i and Sl+i =S\ (periodic 
boundary conditions). In the limit L^oo, the spectrum of this Hamiltonian has been 
exactly solved for at P — 1 using Bethe-ansatz techniques, and the resulting energy 
spectrum is gapless [19, 20, 21]. The ground state of the Hamiltonian has also been 
found [22, 23] for even values of L at P = — 1/3, and the gap above it has been proven to 
be nonvanishing. It has been argued [24] that P — 1 marks a phase transition between 
the gapped Haldane phase [25, 26] and the gapped dimerized phase, and that this phase 
transition is in the same universality class as the two-dimensional Ising model. These 
arguments have been supported numerically [27, 28]. 

Figure 7 shows the estimates S L (P) obtained by exact diagonalization of small 
systems. The finite-size effects are quite strong, and we were unable to make an 
extrapolation of the form (10), since we require < P{l}, c < 1 for all extrapolations 
in order for X and Y to be real and finite. Instead we solve 



L 



T~L — \Sj • Si+i — P(Si ■ S i+ i) 2 



(24) 



i=i 



5 L (P) = D(P) + A(P)L~ B ^ 
5 L ,(P) = D(P) + A(P)L'- B ^ 
5 L „(P)=D(P) + A(P)L"- B ^ 



(25) 



f The symbol A 2 is a part of the name of the numerical method and should not be confused with 
cither A(/3) or A {i} (/3). 
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for B(f3) at fixed (3, and then solve 

5v{P) =A LtI/>L „{l3) + C(P)L 
5 L „(P)=A LtI/jL „(P) + C(p)L 



I-ZB(J3) 



rr-ZB(f3) 



(26) 



for A L)L /, L «(/3) using the value of B{(3) found from (25). Here L < V < L" . The 
parameter Z allows us to tune the extrapolation as in the preceding sections. Figure 8 
shows an example of extrapolations which in fact yield the best estimate of v. 

Even with the extrapolation procedure outlined above, we are not able to eliminate 
the curvature from the CAM plot as we did in sections 2 and 4. Instead, we perform a 
fit to the form 



while varying Z to find u, assigning an arbitrary fixed weight equally to all of the CAM 
points. The best estimate for v is the one that minimizes x 2 [figure 9]. The best fit, 
shown in figure 10, is for i/rs 0.987, in good agreement with theoretical predictions. 

6. Conclusion 

In this paper we propose a new and very general method for constructing canonical 
sequences for use in the Coherent Anomaly Method. A critical point (3 = (3* is 
always marked by the vanishing of some quantity A(f3), though approximations of A(f3) 
often remain nonzero for all values of (3. By intentionally underestimating the rate of 
convergence of an initial sequence of approximations {8i((3)}, we form extrapolations 
{A{j}(/9)} that cross zero at some value of (3, which serves as the approximate critical 
point {(3{i},c}- Furthermore, A{j}(/3) can very generally be expected to cross zero linearly 
with (3. Because there are many ways in which the extrapolations can be made from an 
initial sequence {Si}, we have a great deal of freedom to choose an extrapolation which 
shows a clear coherent anomaly. 

We apply this method to the square lattice Ising model, using the exact values of the 
mass gaps for semi-infinite systems of finite width N at temperature T = f3~ l as 5l(/3)- 
We find that the method yields v L 1 for moderate values of L and that lim^oo v L — > 1, 
as should be expected from the exact solution of the two-dimensional Ising model [5]. 
The convergence may be complicated and nonmonotonic, however, depending on which 
sets of L are chosen to form the extrapolations. 

In order to investigate the convergence of vl to v further, we assume that Sl follows 
the scaling equation (11). This allows us to show that v L ~v + O(L^), which is the 
same convergence rate as has been found for finite-size scaling [11]. 

We also apply this method to the cubic-lattice Ising model, using numerical transfer- 
matrix values of the mass gaps for semi-infinite systems with LxL cross sections and 



Y = — + b + (v-l)X 



(27) 



9 



periodic boundary conditions. Even though we are limited to systems with L < 5 and 
have a limited amount of data for L = 5, we are able to estimate v ~ 0.629, which is 
within a few percent of the best current estimates of v [12, 13, 14]. 

Finally, we apply this method to a one-dimensional spin-1 quantum Heisenberg 
antiferromagnet, using numerical exact-diagonalization estimates of the energy gaps for 
systems of width L and periodic boundary conditions. In spite of large finite-size effects, 
we are able to estimate v ps 0.987, in good agreement with theoretical and numerical 
studies indicating u = l [24, 27, 28]. 

There are a few points that need to be emphasized. 

• Although we use data from finite systems of successively larger size, we are 
not performing Finite-Size Scaling [8]. The initial sequence {<5j} could be 
derived from other techniques, such as the Density-Matrix Renormalisation Group 
algorithm [29, 30], in which the system size is not the most important parameter 
affecting the quality of the approximation. 

• Although we use extrapolations, we are not seeking the best extrapolations in the 
sense of extrapolations which are nearest to the thermodynamic limit of the model 
in question. This is because the extrapolation is just one step in the process. 
Instead, we seek a sequence of extrapolations for which the coherent anomaly is 
clear. 

• Taylor expansions have been used previously to form sequences of approximations 
with i/{j} = l from series expansions [31], but in those studies the linear behaviour 
of the CAM plot could not be improved, since the physics of the series expansion 
left no room for change. The extrapolation-CAM method has flexibility to choose 
the "straightest" CAM plot. 

• This method requires good numerical precision. This is clear from figure 6. Care 
should therefore be exercised when applying this method to Monte Carlo data, 
where statistical uncertainty may be significant. 

• Although we here are estimating only u, other critical exponents could be found in 
the same way. For instance, instead of extrapolating the gap, one could extrapolate 
the reciprocal of the specific heat. 

• Some initial sequences {Si(/3)} may already cross zero as f3 is varied. For instance, 
the variational method of references [32, 33] produces gap estimates for the spin-1 
Heisenberg chain that have this property. In such a case, extrapolations could still 
be used to look for a clearer coherent anomaly. The extrapolations should then 
be faster than the convergence of {Si(P)}, rather than slower, so that the sequence 
{A{j}} still crosses zero. 
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Figure 1. The mass gap, or inverse correlation length, for the square-lattice Ising 
ferromagnet. The dashed lines are the gaps for semi-infinite systems of width 4, 9, 16, 
and 25, from top to bottom. The solid lines are extrapolations using (10) with the 
value of B determined in figure 3; from left to right, they represent {L\,Li2) = (4,9), 
(9,16), and (16,25). 
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Figure 2. CAM plot for the square-lattice Ising model. Equations (8) and (7) define 
X and Y. Extrapolations use the value of B determined in figure 3 and system sizes 
that are consecutive perfect squares: L 6 {4, 9, 16, ... , 400}. Note the difference in the 
scales. 




Figure 3. By varying B, v can be estimated from consecutive extrapolations. In 
order to determine a good estimate of v, we vary B until two consecutive estimates of 
v are identical. This yields 5 = 0.407404 833 and i/ = 0.987 405 623. 
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Figure 4. Using consecutive perfect squares, i.e. Lj = (k + j — l) 2 , j e {1,2, 3, 4} and 
k 6 {1,2,3, . . .}, we determine v and B as in figure 3. As the system sizes become 
large, it appears that j/— >1 and B—>1. For some other series {Lj}, the convergence is 
much less clear. 
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Figure 5. The mass gap for the cubic-lattice Ising ferromagnet. The thin solid curves 
are clamped splines to the transfer-matrix calculations (open circles) of gaps for semi- 
infinite systems of cross section L = 2, 3, 4, and 5, from top to bottom. The thick 
solid curves are extrapolations using (10); from left to right, they represent (Li, L-i) = 
(2,3), (3,4), and (4,5). 
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Figure 6. As in figure 3, we vary B to obtain a good estimate of v. Results are shown 
using two different cubic splines through the data for L = 5. These curves become 
meaningless as (3 C approaches the end of the data available for L — 5. By taking five 
points (solid circles) where the two curves agree and extrapolating, we find v — 0.6285 
(open circle). This compares well with other estimates, such as that of [12] (cross). 




Figure 7. Gap estimates Sl for the spin-1 Heisenberg chain with periodic boundary 
conditions. Strong finite-size effects obscure the fact that the gap vanishes at (3=1 in 
the limit L^oo. 




Figure 9. Wc fit the CAM data to equation (27) while varying Z in (26). The CAM 
data are equally weighted with an arbitrary fixed weight. The minimum of x 2 gives 
the best estimate for v. 




Figure 10. CAM plot for the spin-1 Hciscnberg chain. Y and X are give by (7) and 
(8), respectively. The solid curve is a fit to (27). The dashed line is the asymptote of 
the fitted curve. 



